// Copyright (c) 2009-2010, Jeremy Brewer
// All rights reserved.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are met:
//
//    * Redistributions of source code must retain the above copyright notice,
//      this list of conditions and the following disclaimer.
//    * Redistributions in binary form must reproduce the above copyright
//      notice, this list of conditions and the following disclaimer in the
//      documentation and/or other materials provided with the distribution.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
// POSSIBILITY OF SUCH DAMAGE.

package galaxie500.numeric;

import galaxie500.util.Util;

/**
 * Class which implements a method for adaptively sampling a function
 * 
 * Some functions are expensive to compute (e.g. functions that involve
 * integrals of other functions).  In those cases, we wish to compute the
 * expensive function on a grid of points and interpolate to speed up the
 * computation.  This class implements a Function which does exactly this
 * internally; it uses adaptive sampling to reduce the number of function
 * computations that are necessary.
 * 
 * In order to sample a function, one must define some notion of sampling
 * accuracy.  This class takes its definition of accuracy as being the accuracy
 * in the integral of the function over a given integral.
 * 
 * @author Jeremy Brewer
 */
public class AdaptivelySampledFunction implements Function {

  /** Value of integral */
  private double integral = 0.0;

  /** Internal lookup table */
  private Interpolator interpolator = null;
  
  /**
   * Adaptively samples f over the interval [a, b] to a precision of 1.0e-6.
   *
   * @param f The function to sample
   * @param a The starting point of sampling
   * @param b The end point of sampling
   */
  public AdaptivelySampledFunction(Function f, double a, double b) {
    this(f, a, b, 1.0e-6);
  }
  
  /**
   * Adaptively samples f over the interval [a, b] to the specified precision
   * and builds a lookup table for interpolating the function values.  This
   * constructor is expensive as all of the integration and interpolation is
   * done within it.
   * 
   * @param f The function to sample
   * @param a The starting point of sampling
   * @param b The end point of sampling
   * @param maxError The maximum error to allow in the integral over [a, b]
   */
  public AdaptivelySampledFunction(Function f, double a, double b,
                                   double maxError) {
    ValueStoringFunction g = new ValueStoringFunction(f);

    Integrator integrator = new Integrator();
    integrator.setMaxError(maxError);
    integral = integrator.integrate(g, a, b);

    assert g.xValues().size() > 0;

    // Copy to arrays.
    double[] xValues = Util.sortAndRemoveDuplicates(g.xValues().toArray());
    double[] fValues = Util.sortAndRemoveDuplicates(g.fValues().toArray());

    interpolator = new Interpolator(xValues, fValues);
  }
  
  /**
   * Returns the function value at x if it lies within the original interval
   * passed to the constructor or throws a ValueError if x is outside the
   * interval.
   * 
   * @throws ValueError If x lies outside the interpolated region
   * @param x The value the evaluate at
   * @return An interpolated estimate of the function at x
   */
  @Override
  public double valueAt(double x) {
    return interpolator.valueAt(x);
  }

  /**
   * Returns the value obtained for the initial integration.
   * 
   * @return The integral value.
   */
  public double integral() {
    return integral;
  } 
}
